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M.  Mazumdar 

Westinghouse  Electric  Corporation 

I .  Introduction  to  the  Problem 

In  certain  applied  areas,  notably  that  of  systems  reliability, 
but  also  those  of  work  sampling  and  communications,  it  is  common  to 
consider  a  system  that  randomly  occupies  one  of  two  states.  Thus  in  the 
reliability  application  (upon  which  we  shall  concentrate  hereafter)  a 
radar  system  may  be  operative  ("up"),  or  inoperative  ("down"),  while  in 
work  sampling  a  worker  may  be  "working"  or  "resting",  and  in  communications 
a  telephone  line  or  other  facility  may  be  "idle"  or  "busy".  Random  state 
occupancy  means  that  the  durations  of  the  times  spent  in  the  two  states 
exhibits  random  variability.  For  example,  the  tendency  of  system  times- 
to- failure  (up-times)  to  be  approximately  exponentially  distributed  is  well- 
known  (see  Davis  [  5  ]).  Similarly,  the  idle  time  of  a  telephone  line  will 
be  exponential  if  calls  occur  in  a  stationary  Poisson  fashion,  as  may 
frequently  be  a  reasonable  assumption.  Of  course  distributions  other  than 
the  exponential  can,  and  frequently  do,  arise. 

Our  purpose  in  this  paper  is  to  discuss  the  problem  of  estimating 
the  parameters  of  the  probability  distributions  specifying  the  two-state 
process  described.  Of  interest  also  is  the  estimation  of  functions  of 
these  parameters,  e.g.  "system  readiness"  and  "system  reliability"  in  the 
reliability  context;  these  notions  will  be  defined  later.  The  estimation 
problem  clearly  depends  upon  the  manner  in  which  observations  are  taken, 


2. 


i.e.  upon  the  sampling  plan.  We  assume  here  that  at  least  certain  of  the 
observations  are  available  only  at  rather  infrequent  times.  When  observations 
are  made  we  assume  that  they  either,  a),  reveal  only  the  state  of  the  system 
that  prevails  at  the  instant  of  observation,  in  which  case  we  call  them 
snapshots,  or  b),  consist  of  continuous  recording  of  the  system  state 
throughout  intervals  of  fixed  or  random  duration,  in  which  case  we  call 
them  patches.  In  this  paper  we  shall  consider  mixtures  of  these  plans,  in 
which  snapshots  and  patches  are  mutually  interspersed. 

The  motivation  for  considering  such  sampling  plans  is  essentially 
that  of  determining  system  parameters  economically  and  without  the  bias  that 
is  possible,  at  least  in  the  reliability  area,  from  the  use  of  operator  log¬ 
book  data.  The  latter  frequently  is  extremely  spotty  and  unreliable,  and  a 
supplement  is  desirable.  Snapshot,  or  even  patch,  observations  approximate 
those  made  when  occasional  readings  on  system  state  are  possible.  For 
example,  after  the  time  of  a  military  action  it  is  likely  to  be  accurately 
recalled  whether  or  not  a  radar  was  up  when  needed,  and  possibly  even  the 
approximate  length  of  time  it  remained  so.  See  Cox  ([  1  ],  p.  87  ff.)  for 
a  discussion  of  a  related  machine  stoppage  question.  In  a  recent  paper, 

[  3  J,  Cox  has  also  considered  an  estimation  problem  very  similar  to  ours. 

In  activity  sampling  or  time  study  snapshot  observations  are  frequently 
taken;  it  has  been  suggested  by  Moder  [10  J  that  the  efficiency  of  such 
studies  can  be  increased  by  also  observing  the  remaining  time  in  state  of 
the  system,  particularly  if  the  state  observed  occurs  relatively  infrequently. 
In  telephony  analogous  sampling  procedures  may  also  be  considered. 
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An  interesting  effect  that  sometimes  accompanies  system  observa¬ 
tions  is  that  the  system  behavior  actually  changes  as  a  consequence.  For 
one  thing,  inspections  may  increase  the  vigilance  and  motivation  of  support 
personnel,  thus  inducing  the  system  to  spend  more  time  in  the  up,  or 
desirable,  state.  On  the  other  hand,  too- frequent  observations  are 
distracting  and  time-consuming,  and  may  breed  resentment.  We  shall  not 
attempt  here  to  provide  statistical  procedures  accounting  for  the  possible 
behavioral  impact  of  inspections  or  audits,  although  work  in  this  direction 
is  in  progress. 
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II.  The  Exponential  Model  and  Maximum  Likelihood 

Before  proceeding  to  estimate  parameters  in  the  general  model 
described  it  is  necessary  to  specify  the  latter  mathematically.  In 
reliability  terminology  we  assume  that  the  system  is  alternately  up  and 
down,  that  the  i-th  up  time  is  U.,  the  i-th  down  time  is  D^,  and  that 
[ll^,  i=l,2,  ....]  and  [D^,  i=l,2,  ...]  are  mutually  independent  sequences 
of  independent  and  identically  distributed  random  variables.  Thus  the 
system  is  described  by  a  two-state  renewal  process  (Cox  [  1  ].  )  .  In 
this  case  if  the  system  is  observed  at  widely  spaced  instants  the  chance 
that  it  is  observed  up  is,  by  general  renewal  theory, 

P{up  in  long  run}  =  g[uJ+E^DJ  (2,1) 

with  the  complementary  probability  applying  to  the  complementary  event. 

It  is  thus  clear  that  by  rare  snapshot  sampling  we  can  only  estimate  the 
ratio  of  the  two  means.  However,  if  actual  time  durations  are  observed,  as 
in  patch  sampling,  it  appears  that  estimates  of  both  E[u]  and  E[d]  are 
available.  We  treat  this  possibility  in  terms  of  the  special  model,  in 
which  up-  and  down-times  are  exponentially  distributed. 

A.  Exponential  Model 

Let  the  up-  and  down-times  have  densities 

fu(x)  =  11  >  °>  x  2  °> 

.-Vv  (2.2) 

f D ( y )  -  e  ,yX,  X  >  0,  y  >  0. 


and 
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where  E[u]  =  p  ^  and  E[d]  =  X  *  . 

This  simple  two-state  Markov  model  occurs  in  reliability  applications,  and 
perhaps  in  other  areas  as  well.  We  shall  base  our  estimation  procedures 
upon  it.  However,  we  shall  also  attempt  to  test  the  robustness  of  the 
estimates  thereby  derived  by  treating  processes  having  up  and  down  times 
actually  described  by  other  distributions  (e.g.  members  of  the  log-normal 
family)  as  if  they  were  from  the  Markov  process  generated  by  distributions 
(2.2).  Results  of  this  test,  conducted  by  experimental  sampling,  are 
described  in  Section  IV  . 

B.  Maximum  Likelihood 

Given  the  two-state  renewal  model,  with  distributions  (2.2)  we  can 
now  write  down  the  likelihood  function  under  the  assumption  of  rare  observation. 
Let  us,  for  illustration,  suppose  that  the  system  is  initially  observed  in  the 
up  state,  and  that  starting  at  that  moment  a  patch  of  complete  up  periods 
(U  realizations)  and  complete  down  periods  (D  realizations)  is  recorded. 
Following  this,  a  long  interval  is  allowed  to  elapse  and  the  system  is  again 
observed}  the  state  observed  recorded  at  the  moment  of  first  observation, 
and  a  subsequent  patch  of  up  and  down  realizations  again  recorded;  the 
process  is  repeated  periodically,  always  with  long  intervals  elapsing  between 
consecutive  patches. 

Owing  to  the  assumption  of  long  delays  between  patches  we  shall 
assume  that  the  probability  that  the  system  is  ir.  an  up  condition  when  a  new 
patch  observation  begins  is 
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the  system  is  down  at  patch  beginning  with  the  complementary  probability. 
Moreover,  the  memoryless  or  Markovian  property  of  the  exponential  distri¬ 
bution  (see  Feller  [  6  ]),  assures  that  the  duration  of  the  first  time  in 
state  at  patch  beginning  (termed  "remaining  life",  or  "forward  recurrence 
time",  see  Cox  [  1  ])  is  exponentially  distributed,  with  the  parameter 
appropriate  to  the  state  observed  at  patch  beginning.  It  is  this  last 
property  that  renders  patch  sampling  of  the  exponential  process  especially 
tractable,  for  in  general  the  forward  recurrence  time  is  not  distributed 
in  the  same  manner  as  is  the  corresponding  time  in  state. 

The  above  considerations  then  lead  to  the  likelihood  function 


l(\,p)  - 


e'^V 


e"Xy+\b 


'x+n'  ’ 


(2.4) 


where 


a 

x  =  total  uptime  observed;  x  =  £  x. , 

i=l  1 

x^  being  an  individual  up  interval  (or  forward  recurrence  time), 

a  =  total  number  of  up  intervals  (including  forward  recurrence 
times) , 

y+  =  total  down  time  observed;  y+  =  y^  , 

i=l 

b  =  total  number  of  down  intervals 

(including  forward  recurrence  times), 

a  =  total  number  of  patches  beginning  with  system  up, 

P  =  total  number  of  patches  beginning  with  system  down. 

Expression  (2.4)  can  be  specialized  to  account  for  a  number  of  alternative 
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sampling  plans*  Two  examples  follow. 

Case  li  A  system’s  up  and  down  history  is  continuously  recorded  through 
k  initial  up  and  down  periods.  Thereafter,  m  rare  snapshot  observations 
are  made,  on  r  of  which  the  system  is  in  an  up  condition.  In  this  case 
a  =  b  =  k 
a  =  r,  P  =  m-r  . 


Case  2t  A  system  is  observed  m  times  at  rare  intervals,  and  each  time 
the  system  state  and  the  remaining  time  in  that  state  are  recorded.  Then 
a  =  a  =  r 
b  -  0  =  m-r. 

Returning  to  the  general  case,  with  likelihood  (2.4),  differentia 
tion  of  the  log-likelihood  yields  the  two  equations 


0  log  L  _  a+b  a+P  _  Q 

0X  X  '  \+H  ‘  y+ 


(2.5) 


and 


i  lo9  L  =  a+£  _  o+£ 

OH  V-  X+p  + 


(2.6) 


Eliminating  the  term  involving  (X+n)  ^  we  obtain 


X"1  .  Z£-±  +  a+fi  ,-1 
a+p  a+b  ^  5 


then  substitution  into  (2.6)  produces  the  quadratic  equation  in  p: 


^2[x+(y+“  x+)]  +  p[y+(a-a)  +  x+(2a+b+p )  ]  -  (a+p)  (a+b)  ■  0  . 


(2.7) 
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The  appropriate  solution  is  the  maximum  likelihood  estimate 
a  [y.(a-o)-x  (2a+b+p)]+  *[y+(a-a)-x+(2a+b+p)]2-t4(y+-x+Ha->#)(a+b)x+ 

“ - - ZTTtt - (2>8) 


By  symmetry ,  the  maximum  likelihood  estimate  of  X  is 


A  [x  (b-0)-y.(2b+a+a)]+  'J [x+(b-0)-y+(2b+a+a)]2+4(x+-  y+)(b+a)(a+b)y+ 

X  * - 2Tx+-y>; -  (2-9) 


A  little  algebra  shows  that  the  quantities  under  the  radicals  are  non¬ 
negative.  For  the  Cases  1  and  2  above  (2.8)  and  (2.9)  become 

[y+(k“r)-x+(3k-hn-r)]+  ^  [y+(k-  t)-x+(3k-hn-r)  ] +  8(y+-  x+)  (k*m-r)kx+ 


A 

V-  = 


2(y+-x+)x+ 


(2.10) 


A  [x+(k+r-m)  ~y+(3ktr )  ]+  V  [x+(k+r-m)-y+(3k*r)  ]2*8(x+-y+)  (k+r)ky+ 

- 2Txry;-)7; -  (2.11) 


x  = 


Case  2t 


A 

P  = 


m 


xf  +  \f*jl 


(2.12) 


A 

X  = 


m 


y+ +  V  v+ 


(2.13) 


The  asymptotic  variance-covariance  matrix  can  be  computed  directly. 
We  find  first 


.E[!—L22J=]  .  k-  +  — aA- 

ep  p  p(x+  p) 


(2.14) 


(2.15) 


By  similar  manipulation  there  results  for 


Case_2» 

Var[p] 


X2(Xf2u) 

2mp 


Cov[\,p]  =  ^ 

-1  A 

Under  certain  conditions,  e.g.  when  k  becomes  large,  X  and  p  can  be 
expected  to  be  approximately  normally  distributed}  also,  the  covariance 
tends  to  zero.  A  proof  of  asymptotic  normality  following  the  pattern 
of  Cramer  ([4],  p.  366)  could  be  given,  but  is  omitted. 

Although  the  large-sample  properties  of  maximum  likelihood 
estimators  are  familiar,  it  is  interesting  to  investigate  samples  of 


realistic  size*  We  shall  carry  out  such  an  investigation  here  by 
experimental  sampling.  By  this  means  we  are  able  to  get  an  idea  of 
the  adequacy  of  the  maximum  likelihood  estimators  and  to  evaluate 
simple  empirical  adjustments  to  the  latter  to  improve  their  performance. 
Section  111  contains  such  results.  Sampling  experiments  also  can  be 
expected  to  reveal  possible  inadequacies  of  assumptions  such  as  those 
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made  in  writing  down  the  basic  likelihood  (2.4). 

C.  First  Sampling  Experiment 

This  experiment  was  conducted  for  the  sampling  plan  of  Case  1 
above,  with  estimates  (2.10)  and  (2.11).  The  actual  system  sampled 
involved  parameter  values 

p  *  0.2  (expected  up-time  of  5  units) 

X  ■  1  (expected  down-time  of  1  unit) 
k  *  5  (number  of  initial  periods  observed) 
m  =  10  (number  of  later  snapshots) 

A  synthetic  system  realization  was  observed  continuously  through  five 
consecutive  up-  and  down-times,  after  which  merely  the  state  --  up,  or 
down  —  was  noted;  snapshots  were  taken  at  intervals  of  approximately  15 
time  units.  Five  hundred  realizations  were  examined.  From  the  data  for 
each  realization  estimates  p  and  X  were  computed,  using  (2.10)  and  (2.11). 

In  addition  maximum  likelihood  estimates  of  the  above  parameters  were 
obtained,  using  merely  the  outcomes  of  the  k  =  5  initial  period  observations, 
omitting  the  snapshots;  we  denote  the  latter  by  p  and  X.  Then,  using  both 
sets  of  estimates,  i.e.  (p,  X)  and  (p,  X)  we  computed  estimates  of  system 
performance, 


a)  Operational  Readiness  =  R  =  ^ 

(2.20) 

b)  Operational  Reliability  =  r(T)  -  e  ^ 

(2.21) 

by  substituting  the  estimates  in  for  the  unknown  true  parameters.  The 
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results  estimate  a)  the  long-run  probability  that  the  system  will  be  up 
when  needed,  and  b)  the  long-run  probability  that  the  system  will  be  up 
and  remain  so  for  T  time  units  thereafter,  (e.g.  throughout  a  mission  time). 
For  our  sampling  experiment  T  =  2  time  units. 

A  summary  of  the  results  obtained  is  given  in  the  table  below. 
These  summaries  have  been  computed  as  follows.  Consider  the  estimate  of  p, 
say,  obtained  on  the  i-th  run;  call  it  p^.  Then 


and 


while 


(2.22) 


(2.23) 


(2.24) 


A 

exactly  the  same  procedure  applies  to  \  . 

The  numbers  in  parentheses  next  to  the  estimated  variances  for 
Patch-Snapshot  sampling  were  computed  using  the  asymptotic  formulas 
(2.17)  and  (2.18).  Those  in  square  brackets  next  to  the  corresponding 
means,  variances,  and  mean-square  error  for  Patch  sampling  were  computed 
using  the  exact  formulas 


(2.25) 


Var[f]  -  (jSj)  ^ 


and 


(2.26) 
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■  (k-lllM)  I**  *  (2‘27> 

By  and  large,  agreement  between  the  estimated  variances  and  the  exact 
variances  for  Patch  Sampling  is  very  good.  The  agreement  between  the 
asymptotic  variances  and  sampling  variances  for  Patch-Snapshot  Sampling 
is  less  good.  The  difference  is  probably  attributable  to  the  failure  of 
the  asymptotic  formula  for  our  small  sample  size  (k  =  5,  m  =  10).  An 
alternative  approach  would  be  to  use  the  linearization  or  "delta  method", 
expanding  p  and  X  in  powers  of  x,  y,  and  j-j  and  taking  expectations.  By 
this  means  the  approximate  bias  as  well  as  the  variance  can  be  determined. 


Table  1 

Estimates  of  Failure  and  Repair  Rates 
\=1  ,p=0.2,  k=5,  m=10 

(Statistics  based  on  500  independent  realizations) 
Patch-Snapshot  Patch 


Average si 

A(u) 

=  0.231 

A(p)  - 

0.244 [.250 ] 

a(x) 

=  1.20 

A  (X)  . 

=  1.24  [1.25] 

Variances: 

V(p) 

=  0.011(0.0066) 

V(p)  > 

=  0.020  [.021] 

v(x) 

=  0.41  (0.16) 

V(C)  : 

=  0.61  [0.52] 

Mean  Square: 

M(p) 

=  0.012 

M(p) 

=  0.022  [0.023] 

Error 

m(x) 

}  0.45 

M(X)  : 

=  0.67  [0.58] 

Indications  are  present  in  this  table  that  the  estimates  are  biassed  upwards, 
but  less  so  when  Patch-Snapshot  information  is  available  than  with  Patch 
information  alone;  the  bias  of  the  latter  can  be  computed,  see  (2.25).  The 
difference  between  A(p)  and  p  is  significant  at  the  two-sided  5%  level, 

A 

as  is  that  between  A(X)  and  X.  In  the  following  section  we  discuss  ways  of 
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reducing  the  bias.  In  spite  of  the  presence  of  bias,  Patch-Snapshot 
Sampling  does  considerably  reduce  both  the  variances  and  mean-squared  errors 
of  the  estimates  of  X  and  p  below  those  of  the  estimates  obtained  from 
Patch  Sampling  alone. 

It  is  of  interest  to  study  the  results  of  utilizing  our  estimates 
of  X  and  p  to  estimate  the  measures  R  (see  (2.20))  and  r(T)  (see  (2.21))  of 
operational  performance.  The  following  table  contrasts  the  results  of  using 
Patch-Snapshot  information,  and  Patch  information  alone. 

Table  2 

Estimates  of  Operational  Performance 
X=l,  p=0.2,  k=5,  m=10 

Patch-Snapshot  Patch 


Averages s 

A($)  =  0.822  (0.833) 

A CR)  *  0.815  (0.833) 

A(r)  =  0.533  (0.558) 

A(r)  =  0.526  (0.558) 

Variances: 

V(f!)  =  0.0057  (0.005) 

V(R)  =  0.0094  (.00774) 

V(r)  =  0.016  (0.016) 

V(r)  =  0.022  (.0218) 

Mean-Square 

M(R)  =  0.0058 

M(R)  =  0.0097 

Error 

M(r)  =  0.017 

M(r)  =  0.0226 

A  a 

The  asymptotic  variances  of  R  and  r,  listed  in  parentheses,  were 
obtained  by  delta  methods  for  example,  from  the  approximation 

R  =  R  +  (X-X)  +  (p-p)  , 

there  results,  if  X  =  E[X]  and  p  =  E[p]  , 


(2.28) 


Var  [R]  S  (f  )2var[X]  +  ^  Var[£]  +  2  |S  |B  Cov[X,J]  . 

A  A 

In  fact  we  do  not  have  E[X]  =X  and  E[p]  =  p  for  our  sampling  experiment. 

The  best  estimates  of  E[X]  and  E[p]  available  are  A(X)  and  A(p)  respectively, 

and  expanding  around  these  values  can  be  expected  to  improve  the  approxima- 

a  a 

tions  (although  the  variances  for  R  and  r  in  Table  2  are  in  excellent 
agreement).  The  above  statements  hold  true  also  for  the  estimates  based  upon 
C  and  p  .  When  the  derivatives  in  (2.23)  are  evaluated  and  the  asymptotic 
variances  and  covariance  (see  (2.17),  (2.18),  and  (2.19))  are  substituted 
into  (2.18)  the  results  are 
Case  1» 

Var^R]  *- - - 5 -  (2.29) 

(X+p)  [2mXp+k(X+p)  ] 


and 


Var[r]  =  e  2*lT 


x2  2 
X  u 

(x+p)2 


^2k4-2kT(X-l-|jt)-fT2rmX|^k(X3-u)2l 


2kmXp  +  k2(X+p)2 


} 


(2.30) 


By  setting  in  =  0  in  the  latter  expressions,  variances  for  Patch  sampling 
result.  On  the  basis  of  our  knowledge  of  the  exact  moments  of  X  and  p 
(see  (2.25)  and  (2.26))  a  somewhat  improved  approximation  to  Var[R]  is 
available.  It  is  now  possible  to  expand  around  the  finite  sample 
expectation,  and  to  use  exact  variances.  The  result  is  — —  times 
(2.29),  the  latter  with  m  =  0.  When  k=0  the  sampling  plan  reduces  to 

consideration  of  m  snapshot  observations  only,  and  the  maximum  likelihood 

r  R(1-R) 

estimator  of  R  is  simply  -  ,  with  mean  R  and  variance  ■  * - L .  It  is 

r  7  m  m 

interesting  that  when  k  =  0  in  our  approximate  formula  (2.29)  for 


Var[R]  the  latter  delivers  the  exact  variance  just  quoted. 

Although  Case  2  will  not  be  explored  to  any  great  extent,  the 
form  of  the  estimate  of  R  implied  by  (2.12)  and  (2.13)  is  easy  to  obtain 
and  appears  rather  striking! 


Case  2 


+  JyI 


and 


(2.32) 


Addition  of  the  excess  life  information  thus  reduces  the  variance  of 
our  estimate  of  R  to  exactly  one-half  of  that  of  the  corresponding 
snapshot  estimate.  If  k  Patch  observations  are  alone  available,  then 
(2.29)  provides  that 


Var[R] 


A-  r2(i-r)2 


(2.33) 


The  variances  associated  with  the  above  estimates  should  be  of  aid  in 
selecting  a  sampling  procedure,  if  such  an  option  is  available. 
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III .  Modifications  of  Maximum  Likelihood  Estimates 

Maximum  likelihood  estimates  possess  many  desirable  large-sample 
properties.  If,  however,  one  is  interested  in  estimates  that  minimize 
such  a  plausible  measure  of  loss  or  estimate  ineffectiveness  as  the  mean- 
squared  error,  then  it  is  known  that  modifications  in  the  "raw"  maximum 
likelihood  estimates  are  sometimes  effective.  An  example  of  a  modification 
(to  remove  bias)  is  furnished  by  the  familiar  practice  of  dividing  by 
(sample  size-1)  rather  than  (sample  size)  when  estimating  the  Gaussian 
variance. 

In  this  section  we  use  a  sampling  experiment  to  evaluate  the 
effect  of  some  heuristic  or  empirical  modifications  in  the  Patch-Snapshot 
estimates.  The  modifications  are  suggested  by  first  considering  the  simple 
components  of  the  problem  —  Patch  information  or  Snapshot  information  alone 
and  adapting  the  results  to  the  problem  at  hand. 

A.  Patch  estimate  unbiassing 

The  maximum  likelihood  estimate  of  p,  based  on  a  Patch  sample  of  k, 
~>  k 

+ 

is  biassed.  Since  the  density  of  u+  is  gamma  we  find  directly  that 


Thus  we  can  easily  remove  the  bias  of  p,  and  likewise  of  X,  by  simply 

fj 

multiplying  p  and  X  by  .  This  suggests  the  following  heuristics  for 
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improving  Patch-Snapshot  estimates  (2.10)  and  (2.11) i 

A-l«  Replace  k  by  k-1  in  the  formulas  (2.10)  and  (2.11)  j 

A-2t  Replace  x+  by  x+,  and  y+  by  y+  in  (2.10)  and  (2.11). 


B.  Patch  estimates  with  smallest  mean  square  error 

Consider  estimates  of  the  form  cp,  where  p.  is  the  Patch  Sample 
estimate.  Then,  again  using  the  gamma  distribution  of  x+,  it  may  be  shown 


that 


k-2  ^ 
—  *  = 


k-2 
x , 


is  the  estimate  of  minimum  mean-squared  error.  This  suggests  considering 
the  procedure 

B-li  Replace  k  by  k-2  in  (2.10)  and  (2.11). 

)c 

B-2t  Replace  x+  by  x+,  and  y+  by  Y+  (2.10)  and  (2.11). 


The  above  adjustments  are  but  a  few  of  the  many  possibilities, 
and  are  simply  proposed  for  experimental  investigation.  Their  apparent 
effects  upon  bias  and  mean-squared  error  will  shortly  be  examined.  Before 
doing  this,  we  consider  the  effect  of  modifying  the  part  of  the  likelihood 
function  effected  by  Snapshots  alone. 


C.  Snapshot  estimate  that  minimizes  the  maximum  mean-squared  error 

Suppose  Patch  information  is  ignored.  Then  it  is  known,  see 
Lehmann  [  9  ]»  that  the  estimate  of  a  binomial  parameter  p(=  yyy)  that 
minimizes  the  maximum  mean-squared  error  is,  in  present  notation,  given  by 


m.m. 


r  \fm  _ 1_ 


1+  \Jm  2(1+  \[m) 


(3.2) 
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This  suggests  that  we  replace  r  by 

+  ±  — a—  ,  (3.3) 

1+Vm  2  (l+^m) 

in  the  estimating  formulas  (2.10)  and  (2.11).  However,  the  value  of  this 
modification  can  perhaps  be  questioned  in  advance  if  there  is  evidence  that 
p  is  much  above  one-half j  see  the  comparison  of  risks  of  the  minimax  and  the 
minimum  variance  unbiassed  estimate,  ^  ,  (Lehmann  [  9  ],  p.  4-24).  The 
sampling  experiment  tends  to  bear  out  our  suspicions  for  the  present  example. 

Still  further  modifications  are  suggested  if  prior  information  is 
formally  incorporated  via  Bayes  '  Theorem. 


D.  Snapshot  generalized  maximum  likelihood  estimate,  beta  prior 

Given  Snapshot  information,  and,  in  addition,  a  beta  prior  distri¬ 
bution,  then  the  modified  likelihood  function  appears  as  follows 

L(r;  b)  =  pr+^(l-p)m"r+ft  5  (3.4) 

where  p  =  ,  and  the  generalized  maximum  likelihood  estimate  p^  is  the 

value  of  p  that  maximizes  L(rj  &)i 


r-f  £ 

Pg  m+  +6 

Let  us  suppose  that  rather  diffuse  information  about  the  value  of  p  = 
is  available  from  previous  experience.  Consider  the  priors  with 

1)  ^  -  4,  6  «  1 

2)  *'=8,6  =  2} 


(3.5) 

X 

X+p 


(3.6) 


both  have  modal  values  at  0.8.  While  1)  is  the  more  diffuse,  spreading 
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p-value  probability  more  uniformly  over  [0,  1]  than  does  2),  neither  appears 
to  strongly  pre-judge  the  issue. 

The  heuristic  modifications  suggested  are  then 
D-li  Replace  r  by  m(j^)  in  (2,10)  anc*  (2.11) 

p  |  n 

D-2*  Replace  r  by  * 

The  estimate  (2.32)  can  be  interpreted  in  a  strictly  Bayesian 
manner  as  the  mean  of  a  posterior  distribution  for  p.  L(r;  b)  is 
proportional  to  that  posterior  when  the  prior  is  of  the  beta  form}  we  get 
(2.32)  precisely  when  the  beta  prior  is  proportional  to  p  ^"^(l-p)6”1. 


E.  Generalized  maximum  likelihood  estimate,  beta  prior 

Examination  of  the  likelihood  function  (2.4)  shows  that  a  beta 

prior  for  may  be  applied  directly  to  it  to  create  a  posterior  density 

A.'h*  , 

K .  .« 

for  X  and  p.  If  the  prior  is  proportional  to  p  (1-p)  ,  then  in  the  general 
likelihood  (2.4)  we  only  need  to  change  a  to  a+^»  and  p  to  p+6  to  obtain  the 
posterior.  In  order  to  obtain  the  values  of  X  and  p  that  are  the  modal 
values  of  the  posterior  (generalized  maximum  likelihood  estimates)  the 
above  changes  may  be  incorporated  into  (2.8)  and  (2.9).  In  particular,  for 
Case  1  the  parameters  of  (2.4)  become 
a  -  b  -  k 


o  *  r+  ,  p  *  m-r+ft  *  m+f+6-(r+f) 
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and  the  generalized  estimates  are  obtained  by  replacing  r  by  r+^ and  m  by 
m+^-rt  in  (2.10)  and  (2.11).  We  shall  examine  the  modifications 
E-li  Replace  r  by  r+4  and  m  by  m+5  in  (2.10)  and  (2.11)} 

E-2i  Replace  r  by  r-+B  and  m  by  m+10  in  (2.10)  and  (2.11). 

These  correspond  to  the  beta  priors  of  (2.33). 

The  above  modifications  are  but  a  few  of  those  that  suggest  them¬ 
selves.  For  example,  we  have  not  studied  the  effect  of  introducing  prior 
probabilities  for  X  and  p.  in  the  Patch  part  of  the  likelihood  function  (2.4). 
Natural  priors  for  this  purpose  are  independent  gamma  densities* 

<tW  -  -k^yr  *  (3-7> 


and 


\|j(\)  =  e"^ 


n  » 


(3.8) 


where  the  prior  expectations  and  variances  of  and  X  are  given  by 

-2, 


E[[l]  =  ,  Var  [p]  =  , 


(3.9) 


E[X]  = 


v 

n 


Var [X] 


£m 


(3.10) 


Densities  (2.34)  and  (2.35)  incorporate  conveniently  with  the  part  of  the 
likelihood  involving  e  ^x+  p,ae  X*5  .  Specification  of  the  priors  can  be 
accomplished  by  selecting  values  for  the  above  expectations  and  variances  and 
solving  for  the  parameters  S,  u,  t],  and  v.  Such  specification  ruay,  as  if 
often  suggested  in  Bayesian  analysis,  reflect  subjective  attitudes  as  well 
as  past  performance  of  similar  equipments.  Given  the  parameters,  then  it  is 
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fcisily  seen  that  the  generalized  maximum  likelihood  estimates  for  Patch- 


Snapshot  sampling  result  by  replacing  x+  by  x++K,  y+  by  y++r)»  a  by  a-ki-1, 
and  t  by  b+v-1  in  the  expressions  (2*10)  and  (2.11).  It  appears  that  a 
natural  bivariate  prior  for  Patch-Snapshot  sampling  is  proportional'  to 


e-[uS+\t)] 


(3.11) 


but  in  the  latter  p  and  \  are  not  independent.  Specification  by  choice  of 
parameters  5,  t),  u,  v,  and  w  will  not  be  attempted. 

The  five-hundred  realization  sampling  experiment,  described 
in  Section  II,  C,  was  utilized  to  evaluate  the  effects  of  the  various 
modifications.  The  following  table  summarizes  the  results. 


Table  3 

Sampling  Study  of  Maximum  Likelihood  Estimate  Modifications 

Exponential  Distributions 

Patch-Snapshot  Patch 
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500 


as  an  estimate  of  the 


DISCUSSION  OF  TABLE  3. 

A.  Modification  of  the  rate  estimates 

1.  Consider  first  the  bias*  Use  of 
standard  deviation  of  A ( p.)  together  with  the  normal  approximation  provides 
strong  evidence  that 

A  A 

a)  The  unmodified  Patch-Snapshot  estimate  (and  also  X)  is 
biassed  upwards.  Likewise,  so  are  the  modified  estimates  D,. 

D-l,  D-2,  E-l,  and  E-2.  The  magnitude  of  this  upward  bias  is 

very  much  the  same  for  modifications  other  than  C;  the  bias  of  C 

a  A 

is  relatively  large  for  the  estimate  of  p,  and  small  for  that  of  X. 

b)  The  modified  estimates  A-l,  A-2,  and  in  particular  B-l,  and 
B-2,  are  biassed  downwards.  Modifications  A-l  and  A-2  are 
closer,  on  the  average,  to  the  true  mean  than  the  unmodified 
estimate,  and  also  closer  than  modifications  B-l  and  B-2. 

c)  The  unmodified  Patch  estimate  p,  (and  also  X)  is  biassed 
upwards  to  a  greater  degree  than  the  Patch-Snapshot  estimates. 

As  anticipated,  modifications  A-l  and  A-2  effectively  remove 
this  bias.  Applied  to  Patch-Snapshot  they  tend  to  over-remove 
it  (see  b)),  which  is  not  surprising. 


2.  Consider  the  mean-squared  error. 

A  A 

a)  The  mean-squared  errors,  M(^)  and  M(X),  of  the  Patch-Snapshot 

A  A 

estimates  X  and  jx  seem  to  be  reduced  by  a  factor  of  approximately 
1.65  by  application  of  modifications  A-l  and  A-2.  Modifications 
B-l  and  B-2  have  about  the  same  effect  upon  mean-squared  error, 
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but  they  have  much  greater  effect  upon  variance:  the  B-modifica- 
tions  reduce  the  variance  of  the  estimates  by  approximately  a 
factor  of  3. 

b)  Modification  C  has  a  more  beneficial  effect  upon  mean  and 

A  A 

mean-squared  error  of  X  than  upon  p.  Modifications  D-l  and 

D-2  have  only  very  slightly  beneficial  effects,  again  more 

A  /\ 

importantly  upon  X  than  upon  p.  The  same  is  true  of  E-l  and  E-2. 

c)  Of  the  modifications  considered,  the  B-modif ications 
produce  the  smallest  mean-squared  errors  for  the  Patch  estimates. 

The  A-modifications  are,  however,  effective  in  reducing  mean- 
squared  error  and  are  exactly  unbiassed. 

B.  The  effect  of  the  modifications  on  estimates  of  probabilities  R  and  r 
1.  Consider  the  bias. 

a)  Use  of  the  sample  variance  leads  to  the  conclusion  that  R  and  R, 
the  Patch  Snapshot  and  Patch  estimates  of  R,  are  biassed  downwards. 
None  of  the  suggested  modifications  significantly  effect  this  bias, 
save  for  C.  which  somewhat  increases  its  magnitude.  Very  probably 
the  fact  that  the  priors  involved  in  the  D  and  E  modifications  have 
a  mode  below  0.833  (at  0.8)  has  a  tendency  to  further  increase 

A  w 

the  bias.  The  bias  of  R  is  somewhat  smaller  than  that  of  R. 

b)  The  Patch-Snapshot  estimate  r  appears  to  be  biassed  downwards 
when  no  modification  is  made,  but  use  of  the  sample  variance  does 
not  enable  us  to  state  this  with  high  confidence.  Again  r  has  a 
greater  bias  than  r.  The  A  modifications  seem  to  remove  this  effect, 
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perhaps  replacing  it  by  a  smaller,  upward,  bias.  The  B  modifi¬ 
cations  apparently  go  too  far. 

2.  Consider  the  mean-squared  error. 

a)  Only  the  D  and  E  modifications  seem  to  have  a  beneficial  effect 

A 

upon  the  mean-squared  error  of  R}  these  effectively  cut  the  mean 
squared  error  in  half.  A  smaller,  but  still  desirable,  impact  upon 
is  noticeable. 

b)  If  sampling  were  Snapshot  alone,  then  the  estimate  of  R  is  the 

r  — 

number  of  successes  —  =  R  ,  and 

m  * 

£[R]  '  &  ’  Var[r]'W' 

Thus  R  is  unbiassed,  and  for  our  example 

V(R)  *  M(R)  =  0.0139  .  (3.12) 

A  <-/ 

Comparison  with  V(R)  and  V(R)  indicates  that  Patch  Snapshot  reduces  the 


above  mean  square  error  by  about  a  factor  of  two. 
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IV.  Robustness 

The  exponential  assumptions  (2.2)  made  in  deriving  our  maximum 
likelihood  estimates  are  often  plausible,  especially  in  reliability 
problems.  However,  one  seldom  can  test  such  assumption:  thoroughly,  and  it 
seems  worthwhile  to  investigate  the  behavior  of  our  estimates  when,  in  fact, 
other  —  perhaps  equally  plausible  distributions  govern  the  observations. 
To  this  end  we  have  utilized  a  sampling  experiment,  first  generating  five 
hundred  realizations  for  each  of  three  sets  of  alternative  distributions, 
then  sampling  from  these  realizations  as  was  done  before,  and  finally 
computing  estimates  of  p,  X,  R,  and  r.  Modifications  were  also  considered. 

Specifically,  suppose  that  the  up  times  of  the  two-state  process 
have  the  log-normal  distribution! 


P  U  <  x 


(In  x-m  )/o 
u  '  u 


le'  Z/*'  dz 


(4.1) 


V2"  i 

and  that  the  down  times  also  have  the  log-normal  distribution  with  para- 
2 

meters  m^  and  .  We  can  express  the  moments  of  U  (and  D)  in  terms  of 
2  2 

m  (m.)  and  o  (o,)}  and  vice-versa,  thus 
u  dud  ’ 


m 


u 


(4.2) 


2 

Given  ■‘•[U]  and  E[U  ],  and  the  corresponding  moments  for  down  time  one  can 
utilize  (4.1)  to  draw  appropriately  centered  and  scaled  random  normal 
deviates.  Let  U  be  such  a  number}  then  a  realization  of  U  having  desired 
mean  and  variance  is  achieved  by  computing 
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From  such  log-normally  distributed  realizations,  a  realization  of  the 
up-down  process  history  may  be  obtained. 

As  was  remarked  earlier,  one  use  for  such  histories  is  to 
investigate  the  behavior  of  maximum  likelihood  estimators  computed  under 
the  exponential  specification  when,  in  fact,  the  log-normal  specification 
prevails.  The  results  of  such  investigation  are  tabled  and  discussed 
below.  We  take,  as  in  the  exponential  situation, 


E  [U]  =5 


E[D]  =  1. 


(4.3) 


Now  for  the  exponential  distribution  the  following  relationship  between 
first  and  second  moments  exists: 

£[U2]  =  2  E2[U]  .  (4.4) 

We  can  easily  determine  log-normal  distributions  with  this 
property.  From  (4.2)  the  underlying  normal  must  in  this  case  have 
parameters 


mu  =  log  E [U ]  -  ~  log  2,  o2  =  log  2.  (4.5) 

A  log-normal  distribution  with  parameters  (4.5)  will  agree  with  an 
exponential  with  mean  E[U]  up  to  second  moments.  In  what  follows  we  shall 
refer  to  this  situation  as  Case  E,  (E  =  exponential). 

Contrasted  with  the  above  will  be  examined  one  in  which  the 
coefficient  of  variation  is  one-half: 

E[U2]  =  |  E2[U], 


(4.6) 
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so 

mu  -  log  E[U]  -  ^  log  (|)  ,  o2  =  log  \  .  (4.7) 

This  is  called  Case  H  (H  *  half-exponential). 

The  last  case  to  be  considered  is  that  in  which  the  coefficient 
of  variation  is  twice  that  of  the  exponential* 

E[U2]  =  3  E2[U]  ,  (4.8) 

so 

mu  «  log  E[U]  -  ^  log  3,  o2  =  log  3.  (4.9) 

This  is  called  Case  T  (T  =  twice-exponential ) . 

Comparisons  of  the  behavior  of  our  estimates  when  the  distribu¬ 
tions  of  Cases,  E,  H,  and  T  prevail  are  given  in  the  next  tables. 


Table  4 

Sampling  Study  of  Maximum  Likelihood  Estimates  of  Rates  X.  and 

Modifications  and  Robustness 
Exponential  and  Log-Normal  Distribution 
Patch-Snapshot  Patch 
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DISCUSSION  OF  TABLE  4. 


a)  Use  of  V  500^  and  the  normal  approximation  provide  strong 

evidence  that  for  Patch-Snapshot  sampling  the  L.-N.  cases  E,  T,  and  H  are 
all  biassed  upwards,  just  as  is  the  original  exponential.  The  upward  bias 
seems  actually  to  be  somewhat  smaller  for  the  cases  L.-N.,  E  and  L.-N., 

H  than  for  the  exponential  case  itself.  If  the  parameters  were  estimated 
knowing  that  the  L.-N.  specification  were  correct  still  better  properties 
would  be  expected  for  the  estimates,  but  there  are  indications  that  only 
when  L.-N.,  T  prevails  does  any  degradation  occur  in  mean-squared  error. 

The  latter  effect  is  possibly  attributable  to  the  positive  skewness  of 
L.-N.,  T,  which  is  greater  than  that  for  the  exponential.  The  degree  of 
bias  exhibited  by  the  Patch-Snapshot  estimates  seems  also  to  be  smaller 
than  that  of  the  Patch  estimates  alone. 

The  mean-squared  errors  of  our  estimates  (2.10)  and  (2.11)  in  the 

A  * 

L.-N.  cases  seem  satisfactory;  for  only  L.-N.,  T  do  M(p)  and  M(\)  exceed 
the  values  for  the  exponential  case.  Again  an  improvement  over  the  mean- 
squared  errors  for  Patch  estimates  alone  seems  to  exist,  although  this 
improvement  is  not  dramatic. 


b)  For  the  L.-N.  specification  the  A-modif ication  appears  to  over-reduce 
the  upward  bias  mentioned  in  a)  above,  at  least  for  L.-N.,  E  and  H.  It  is 
nearly  correct  for  L.-N.,  T.  The  net  effect  is  to  reduce  or  leave 
essentially  unchanged  the  bias,  except  for  L.-N.,  H;  in  the  latter  case 
the  bias  is  increased.  But  in  all  cases  the  mean-squared  error  is  brought 
down  by  the  A-modif ication,  sometimes  considerably. 
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c)  The  B-modification  tends  to  over-correct  the  bias,  without  producing 
a  compensatory  effect  upon  mean- squared  error.  In  case  L.  H  the  mean- 
squared  error  seems  to  have  been  increased. 

Lastly,  we  examine  the  maximum  likelihood  estimates  of  R  and  r  when 
the  log-normal  is  the  actual  underlying  distribution.  A  summary  appears  in 
Table  5. 


Table  5 

Sampling  Study  of  Maximum  Likelihood  Estimates  of  Probabilities  R  and  rj 

Modifications  and  Robustness 
Exponential  and  Log-Normal  Distributions 
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DISCUSSION  OF  TABLE  5. 

a)  There  is  an  apparent  tendency  for  the  estimates  R  and  r  to  be  biassed 
downwards;  the  bias  of  the  latter  being  greater  than  that  of  the  former. 
Notice  that  the  biases  are  actually  not  as  serious  for  L.-N.  cases  as  for 
the  exponential.  The  various  modifications  of  the  estimates  of  X  and  p 
seem  to  have  very  little  effect  upon  the  magnitude  of  tho  bias.  Modifica¬ 
tions  A,  B,  and  C  tend  to  somewhat  increase  the  mean-squared  error  of  R. 

The  D  and  E  modifications  somewhat  reduce  the  mean-squared  error  of  R. 
Certainly  the  latter  mean-squared  error  is  brought  considerably  below  that 
of  the  "distribution- free"  estimate  (3.12). 

b)  The  bias  of  the  unmodified  estimates  r,  and  likewise  of  r,  is  in  a 
downwards  direction.  In  part  this  may  be  the  effect  of  the  downward  bias 
of  the  R  estimates,  and  in  part  that  of  the  upward  bias  of  the  estimates  of 
[i.  Apparently  the  A  and  B  modifications  over-compensate  for  this  effect. 
For  the  present  numbers  the  simple  arithmetic  average  of  the  unmodified 
and  the  A-modified  estimates  has  smaller  bias  than  either  component  alone; 
this  may  be  an  accident,  but  further  heuristics  are  suggested  as  a 
consequence.  Modifications  D  and  E  force  the  bias,  and  a] so  the  mean- 
squared  error,  still  further  downwards.  It  would  be  of  interest  to  combine 
modifications,  applying  A  and  D  together  for  example,  but  this  experiment 
has  not  yet  been  conducted. 

Tables  4  and  5  suggest  that,  at  least  for  the  particular  situation 
considered  in  our  sampling  experiment,  the  maximum  likelihood  estimates 
based  on  an  exponential  specification  behave  well  also  when  the  data  is 
log-normal. 
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Certainly  no  hard  and  fast  conclusions  are  possible  as  a  result  of 
the  sampling  experiments.  However,  possibly  useful  indications  concerning 
the  behavior  of  our  maximum  likelihood  estimates  are  obtained  where  none  seem 
otherwise  available. 
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